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Abstract 

In this preliminary report, a multi- region radial integration is compared to 
the recently proposed method due to Handy et al.. Preliminary results for small 
systems indicate that the new integration scheme is generally comparable to 
and sometimes better than that of Handy et al., although this conclusion is by 
no means firm. Work for larger systems is continuing. 
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1 Introduction 


There is recently a great deal of interest in the testing and application of Density 
Functional techniques for use in molecular calculations. The complicated functional 
forms used in this theory preclude analytical integration, so that numerical methods 
must be used. Recently, Handy et a/.[l] proposed a radial integration based on the 
Euler-Maclaurin summation giving good results for a number of systems. However, 
it is worthwhile investigating other methods since for accurate calculations many 
quadrature points are required in general, even with the Euler-Maclaurin technique. 
In addition, on problem with the Euler-Maclaurin scheme is that there is little control 
over the location of additional quadrature points. For example, this scheme is very 
efficient for atoms, but in a molecular context additional points are required for an 
integration which is comparable to the atomic case. These additional points should 
ideally be placed in the valence region. However, with the Euler-Maclaurin scheme 
points must be placed throughout the entire region. This problem is bypassed with 
a multiregion integration. 

The advantage of the multiregion integration is the control of the placement 
of additional quadrature points in the molecular situation. In general, we determine 
an accurate quadrature for the constituent atoms or small molecular fragments. In 
the molecular situation we then only need to add points to the valence region(s), and 
increase the accuracy of the angular quadrature, in order to be able to produce an 
accurate molecular quadrature. One disadvantage of the multiregion integration is 
a possible increase in errors due to the accumulation of errors from the individual 
regions. 

Here we investigate the efficiency of a multi-region integration similar to that 
previously proposed by Te Velde and Baerends [2] and also Thakar [3]. One difference 
in the current work is the use of a change of vaxiables due to McLean and Yoshimine [4] 
which allows for efficient placement of the quadrature points by mapping the region 
[a, m, 6] onto the region [—1,0, 1], where a and b are the endpoints of the region, and 
m is any point between a and b. This mapping takes the form 

X = Ci+ c 2 ( 1 + p)/{ 1 - fit). 

The constants Ci, c 2 and /3 are completely determined by the mapping a — ► — 1, 
m — ► 0, and b — * 1. This mapping may be used for all ranges except the doubly 
infinite range (-oo, oo), providing appropriate limiting forms are used for special cases. 

The regions used here are based on the atomic charge density, in a spirit similar 
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to that used by Thakar [3]. However, rather than using a general rule based on the 
atomic number by which to determine the regions, we use the atomic charge density 
directly. The function r 2 * p(r) is generated numerically, and analysed for maxima, 
minima and points of inflection. In general, the regions we define are based on these 
points, in particular the points a and b map onto two minima in the charge density, 
and the point m maps onto the maximum in between these. These three points define 
a shell, and the mapping ensures an optimal distribution of quadrature points in the 
region. Generally we have one region for each shell, with the valence shell being split 
into an inner and outer parts at the point of inflection of the valence shell (see also 
Thakar), and the innermost (core) region into two at the first point of inflection. 
Although this scheme may seem somewhat complicated, in fact it is straightforward 
to generate the required regions and tabulate them for any atom. 

For the inner regions we use a standard Gauss- Legendre quadrature, whereas 
for the outer region we use a shifted Gauss-Laguerre scheme based on the weight 
function e~ aT , where a is related to the highest occupied MO eigenvalue (or the IP) 
by a = 2(\J—2ti) (see also Thakar [3]). 

The other details of the integration are very similar to that proposed by 
Handy et al. [1]. We use a product scheme for the angular integration identical 
to that of Handy et al., including the use of the angular crowding parameter K$. The 
single center integration scheme of Becke [5] is used, but we use the cutoff function 
due to Handy et al., with standard Bragg-Slater radii. 

We consider three systems here - Ne, Cu and CO. We compare the Euler- 
Maclaurin to the current multiregion integration for all four systems. We consider 
integration of the total charge density only here, at the SCF level of approxima- 
tion, since Handy et al. [1] showed that this gives a good indication of the overall 
performance of an integration scheme. For Ne, we use the Dunning [6] [5s4p] con- 
tracted Gaussian basis set. For carbon and oxygen we use the correlation consistent 
polarized valence triple zeta (cc-pVTZ) basis sets of Dunning. For Cu we use a 
(20s 15p lOd 6/ 4 p)/[6 + Is 5 + lp id 2/ lp] basis set derived from the large primi- 
tive set of Partridge [7], contracted by Bauschlicher using the atomic natural orbital 
(ANO) approach. The C-0 bond length is 2.2 a.u. in all cases. 

The calculations were carried out using the Seward integral program and the 
Sweden SCF program on the Cray YMP-C90 at NASA Ames research center. 
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2 Results and Discussion 


For Ne atom, the results for the Euler-Maclaurin integration are given in Table 1. 
We give the number of points used in the radial integration, the number of points 
including those discarded due to the radial cutoff (in parentheses), the radial factor 
m T (Handy et al. [1], Eqn. 6), and the total error in the integrated charge density. 
Our results for this are very similar to Handy et al., as expected, with a very accurate 
radial integration attained with 72 radial points. In Table 2, we give the results from 
the multiregion integration, with the number of points in each region given. In region 
5 we also indicate the number of points including those discarded due to the radial 
cutoff, in parentheses. Overall, more radial points are required for the multi-region 
integration, although the integration is still accurate. For example, if one used a 
value for m T in the Euler-Maclaurin integration which was not optimal, then the 
results could be worse than for the multi-region integration. 

In Table 3, we give the results for the integration of Cu atom in a large ANO 
basis set, for three different numbers of radial points. One trend to be noted is that a 
higher m T value is needed for Cu than Ne as more radial points are added, and more 
radial points are needed to achieve a similar absolute accuracy. As noted by Handy et 
al. this is because the numerical integration only gives a certain relative error rather 
than absolute error in the integration. We note that for Cu atom (4s 1 ) it is necessary 
to integrate out a long way (22 a.u.) due to the very diffuse nature of the 4s orbital, 
and that sometimes a higher value of m T can be necessary than that recommended 
than Handy et al. to obtain optimum results. 

The multi-region results of Tables 4 and 5 compare quite well with the Euler- 
Maclaurin results, although again it seems that more points are necessary for the 
multi-region scheme, depending on the m T value used for in the Euler-Maclaurin 
scheme. We note that the inner regions for the multi-region scheme are very compact, 
due to the higher atomic number in this case, and that quite a few radial points are 
necessary to describe the density accurately in this compact region. However, this 
does not necessarily translate to a lot of points in a molecular calculation, since the 
number of angular points needed in the inner regions is much smaller than in the 
valence regions. 

Finally we consider the CO molecule (Tables 6-8). We use the same num- 
ber of angular points in for the Euler-Maclaurin and the multi-region integrations 
{n$= 42, n^=84), and the same Bragg-Slater radii in each case (rc=l .32281 a.u., 
ro=l. 13383 a.u.). Overall, it seems that it is possible for the Euler-Maclaurin scheme 
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to outperform the multi-radial scheme if the right values of m T and (the angular 
factor [1]) are chosen. However, these differ significantly from the “standard” values 
recommended by Handy et al., m r = 2 or 3 and m^=10 or 11. For these values, we 
see that the multi-region scheme is either equivalent to or better than the Euler- 
Maclaurin scheme in efficiency. Another interesting point from Tables 6 and 8 is that 
it is not just the total number of radial points which is important, but also the spread 
of these points along r, as this affects the total number of points through the angular 
factors. For example, the spread of radial points in the Euler-Maclaurin scheme is 
approximately linear on a logarithmic scale, apart from the very short and long range 
regions. Changing the factor m r changes the slope of this logarithmic plot, so that 
for higher values of m r there are more points in the core region and more points in 
the long range regions, leading to fewer points overall, since there are fewer angular 
points in the core region and the long range points are discarded due to the radial 
cutoff. Thus for the Euler-Maclaurin scheme, very different total numbers of points 
are realized for m r = 3 versus m T = 4, even though the nominal number of radial points 
is the same in each case. 

For the multi-region scheme, two sets of results are presented for CO (Table 8), 
differing in the number of points in the valence and outer regions, and with the total 
number of points being very comparable to the best results found for the Euler- 
Maclaurin scheme. One advantage of the multi-region scheme is that we are able to 
take a set of integration parameters from a similar atom (for example, those for Ne in 
this case), and then place more radial points in the valence and outer valence regions 
in order to attain higher accuracy. This can be seen to be an effective way to add 
points for CO. 

Overall, it seems that neither integration scheme is clearly superior in the 
molecular situation, based on the current results. More and larger systems need to 
be studied in order to establish whether the multi-region integration scheme proposed 
here is significantly better to the widely used Euler-Maclaurin scheme of Handy et 
al.[l] 
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Table 1: Euler-Maclaurin integration of Ne atom (r max =15.0 a.u.) 


n r 

m r 

Ap{r) 

66(70) 

i 

3.0xl0 -5 

56(70) 

2 

3.1xlO~ 10 

50(70) 

3 

8.8xlO -10 

47(70) 

4 

7.1xl0 -7 

94(100) 

1 

6.6xl0 -6 

72(100) 

2 

9.1xl0 -13 

80(100) 

3 

3.8xl0 -11 

67(100) 

4 

9.1xl0 -11 


Table 2: Multiregion integration of Ne atom (r mor =15.0 a.u.) 




n r (region) 


Total n T 

Ap(r) 

1 

2 

3 

4 

5 



mm 

15 

12 

16 

3(3) 

56(56) 

5xl0 -9 

a 

20 

24 

24 

11(16) 

94(99) 

lxlO -12 
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Table 3: Euler-Maclaurin integration of Cu atom (r max =22.0 


n r 

m T 

Ap(r) 

47(64) 

2 

6-OxlO" 6 

43(64) 

3 

6.4xl0" 8 

40(64) 

4 

7.6xl0 -6 

71(96) 

2 

5.2xl0 -9 

64(96) 

3 

1.4xl0 _1 ° 

60(96) 

4 

4.4xlO -10 

95(128) 

2 

6.5xl0 -9 

86(128) 

3 

4.9xl0 -12 

80(128) 

4 

l.lxlO -13 

77(128) 

5 

l.lxlO -10 


Table 4: Regions for Cu atom (a.u.) 



a 

m 

b 

1 

0.000 

0.005 

0.010 

2 

0.010 

0.038 

0.077 

3 

0.077 

0.170 

0.380 

4 

0.380 

0.630 

0.940 

5 

0.940 

2.500 

6.000 

6 

6.000 

— 

22.000 
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Table 5: Multiregion integration of Cu atom (r max — 22.0 a.u.) 


1 

2 

n r ( region) 

3 4 5 

6 

Total n T 

Ap(r) 

15 

15 

15 

10 

20 

15(9) 

84(90) 

2.2xl0 -9 

20 

15 

15 

10 

20 

11(20) 

91(100) 

6.7xlO -11 

25 

20 

20 

15 

25 

11(20) 

117(130) 

1.4xl0 -12 
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Table 6: Euler-Maclaurin integration of CO molecule (ng=42, 

radbs=l. 32281,1-13383 a.u.n^,=84, r maE =15.0 a.u.). Note: same value of m M and m r 
used for C and 0. 


n T 

Ngrid 

m r 


A^(r) 

132,134 (192,192) 

60902,62666 

3 

10 

1.3xl0 -9 

132,134 (192,192) 

60902,62666 

3 

11 

3.0xl0 _1 ° 

132,134 (192,192) 

60902,62666 

3 

12 

2.8xlO -11 

132,134 (192,192) 

60902,62666 

3 

13 

9.8xl0~ 12 

132,134 (192,192) 

60902,62666 

3 

14 

2.5xl0 -11 

148,150 (192,192) 

87162,88926 

2 

10 

2.1xl0 _1 ° 

148,150 (192,192) 

87162,88926 

2 

11 

8.7xl0 -11 

148,150 (192,192) 

87162,88926 

2 

12 

3.9xl0 -11 

148,150 (192,192) 

87162,88926 

2 

13 

2.2X10* 11 

148,150 (192,192) 

87162,88926 

2 

14 

1.6xl0 -11 

148,150 (192,192) 

87162,88926 

2 

15 

1.7xl0 -11 

124,125 (192,192) 

47228,48110 

4 

11 

8.7xl0 -9 

124,125 (192,192) 

47228,48110 

4 

12 

5.4xl0 -9 

124,125 (192,192) 

47228,48110 

4 

13 

1.8xl0' 9 

124,125 (192,192) 

47228,48110 

4 

14 

5.2xlO" 10 

124,125 (192,192) 

47228,48110 

4 

15 

1.5xl0 -9 
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Table 7: Regions for CO molecule (a.u.) 



a 

m 

b 

1 

0.00 

c 

0.025 

0.05 

2 

0.05 

0.17 

0.63 

3 

0.63 

1.27 

1.90 

4 

1.90 

4.0 

6.0 

5 

6.0 

— 

15.0 

1 

0.0 

0 

0.02 

0.04 

2 

0.04 

0.12 

0.41 

3 

0.41 

0.83 

1.37 

4 

1.37 

4.0 

6.0 

5 

6.0 

— 

15.0 






Table 8: Multiregion integration of CO molecule (r mar — 15.0 a.u.) 




n T (region) 


Total n r 

Nprt'd 


Ap{r) 

1 

2 

3 

4 

5 





15 

20 

24 

24 

7(16) 

90(99) 

58078,56944 

12 

3.0xl0 -9 

15 

20 

24 

24 

7(16) 

90(99) 

58078,56944 

13 

l.OxlO -10 

15 

20 

24 

24 

7(16) 

90(99) 

58078,56944 

14 

1.4xl0 -9 

15 

20 

24 

24 

9(24) 

92(107) 

59842,58708 

13 

8.3xl0 -u 

15 

20 

24 

24 

11(32) 

94(115) 

61606,60472 

13 

8.3xl0 -11 

15 

20 

24 

32 

9(24) 

100(115) 

66898,65764 

13 

5.7xl0 -11 

15 

20 

24 

40 

9(24) 

108(123) 

73954,72820 

13 

l.OxlO -11 

15 

20 

24 

40 

9(24) 

108(123) 

73954,72820 

12 

1.4xlO -11 

15 

20 

24 

40 

9(24) 

108(123) 

73954,72820 

14 

1.2xl0" u 




